function  SlicePlot(handles, position)
if cuda %&& false
    CudaSlicePlot(handles,position);
    drawnow
else
    % This function plots the 6 slices as subplot in the panel pointed by
    % panelHandle
    panelHandle = handles.plotPanel;
    posHandle = handles.cursorPosition;
    stack = get(panelHandle,'userdata');
    if(isempty(stack)) return; end;
    % Create subplots at first call / retrieve handles otherwise
    children = (get(panelHandle,'children'));
    children = flipud(children);
    %posHandle = children(1);
    set(posHandle,'data',position);
    %compute Talairach and MNI coordinates
    M_tal = get(handles.talairachTable,'userdata');
    if ~isempty(M_tal)
        [Pac,Pt] = Talairach(position,M_tal);
        set(handles.talairachTable,'data',Pt);
        set(handles.mniTable,'data',Pac);
    end
    if ~isempty(get(handles.selectExcisionVolume,'tooltipstring'))    
        pos = get(handles.cursorPosition,'data');
        pos = [pos 1];
        stack = get(handles.plotPanel,'userdata');
        M = stack(8).mat;
        pos = M\pos';
        V = stack(8).vol;
        pos = int16(pos(1:3));
        try
            N = V(pos(1),pos(2),pos(3));
            load([pwd '\atlas\AreaList.mat']);
            if N~=0
                name = data{N};
            else
                name = '';
            end
            set(handles.selectedArea,'string',{'Area:';name});
        catch
            set(handles.selectedArea,'string',{'Area:'});
        end
    end

    if isempty(children);
        sp = struct('z',[subplot(2,3,1,'Parent',panelHandle) subplot(2,3,4,'Parent',panelHandle)],...
                    'y',[subplot(2,3,2,'Parent',panelHandle) subplot(2,3,5,'Parent',panelHandle)],...
                    'x',[subplot(2,3,3,'Parent',panelHandle) subplot(2,3,6,'Parent',panelHandle)]);
        for j=1:2
            set(sp.x(j),'position', [0.67 1-0.5*j 0.32 0.5]);
            set(sp.y(j),'position', [0.34 1-0.5*j 0.32 0.5]);
            set(sp.z(j),'position', [0.01 1-0.5*j 0.32 0.5]);
        end
    else
        sp = struct('z',children(1:2),'y',children(3:4),'x',children(5:6));
    end

    for (j=1:2) % plot up first then bottom with only IRM
        if j==2
            N = 1;
        else
            N = length(stack);
        end

        cla(sp.x(j));
        cla(sp.y(j));
        cla(sp.z(j));
        for(i=1:N) % draw every images
            if i==7
                continue;
            end
            pixPos = stack(i).mat \ [position 1]';
            P = spm_imatrix(stack(i).mat);
            scaleFactor = P(7:9);
            volSize = size(stack(i).vol);
            % get corner coordinates in real world units.
            corners = zeros(4,8);
            for k1 = 0:1
                for k2 = 0:1
                    for k3 = 0:1
                        corners(:,k1*4+k2*2+k3+1) = stack(i).mat * [k1 * volSize(1) ; k2 * volSize(2) ; k3 * volSize(3) ; 1];
                    end
                end
            end
            MIN = min(corners');
            MAX = max(corners');
            clearance = MIN;
            % Draw z slice.
            matz =  stack(i).mat \ makehgtform('translate',[clearance(1) clearance(2) position(3)]);
            imSize = MAX - MIN;
            imSize = [imSize(1) imSize(2)];
            interpolationMethod = 1;
            if (i==6 || i== 8 || i==9 || i==10)
                interpolationMethod = 0;
            end
            if i~=8
            im = uint16(spm_slice_vol(stack(i).vol,matz,imSize,interpolationMethod));    
            %apply brightness and contrast
            im = ApplyContrastBrightness(im);
            else
                coordX = MIN(1):MIN(1)+imSize(1);
                coordY = MIN(2):MIN(2)+imSize(2);
                im = uint16(zeros(imSize));
                for k=1:imSize(1)
                    for l=1:imSize(2)
                        [xxx talcoord] = Talairach([coordX(k) coordY(l) position(3)],M_tal);
                        talcoord = (stack(i).mat \ [talcoord 1]');
                        talcoord = int16(talcoord(1:3));
                        talcoord = max(talcoord,1);
                        im(k,l) = stack(i).vol(talcoord(1),talcoord(2),talcoord(3));
                    end
                end
            end
            if j==2 alpha = 1; else alpha = stack(i).alpha; end
            alpha = (im' ~= 0) * alpha;
            RGBim = ind2rgb(im',stack(i).colormap);
            image('parent',sp.z(j),'XData', [MIN(1) MAX(1)],'YData', [MIN(2) MAX(2)],...
                  'CData',RGBim,'AlphaData', alpha, 'ButtonDownFcn',{@SetPosition,'Z',position(3),handles});
            % Draw y slice.
            T = makehgtform('translate', [clearance(1) clearance(3) position(2)]);
            T = makehgtform('translate', [clearance(1) position(2)  clearance(3)]);
            RX = makehgtform('xrotate',pi/2,'zrotate',0);
            maty = stack(i).mat \ T * RX;
            imSize = MAX - MIN;
            imSize = [imSize(1) imSize(3)];
           % im = stack(i).vol(:,:,int32(pixPos(1)));
            im = uint16(spm_slice_vol(stack(i).vol, maty, imSize,interpolationMethod));
            im = ApplyContrastBrightness(im);
            if j==2 alpha = 1; else alpha = stack(i).alpha; end
            alpha = (im' ~= 0) * alpha;
            RGBim = ind2rgb(im',stack(i).colormap);
            image('parent',sp.y(j),'XData', [MIN(1) MAX(1)],'YData', [MIN(3) MAX(3)],...
                  'CData',RGBim,'AlphaData', alpha,'ButtonDownFcn',{@SetPosition,'Y',position(2),handles});
            % Draw x slice
            T = makehgtform('translate', [position(1) clearance(2) clearance(3)]);
            RY = makehgtform('yrotate',-pi/2);
            matx = stack(i).mat \ T * RY;
            imSize = MAX - MIN;
            imSize = [imSize(3) imSize(2)];
            
           % im = stack(i).vol(:,:,int32(pixPos(1)));
            im = uint16(spm_slice_vol(stack(i).vol, matx, imSize,interpolationMethod));
            im = ApplyContrastBrightness(im);
            if j==2 alpha = 1; else alpha = stack(i).alpha; end
            alpha = (im ~= 0) * alpha;
            RGBim = ind2rgb(im,stack(i).colormap);
            image('parent',sp.x(j),'XData', [MIN(2) MAX(2)],'YData', [MIN(3) MAX(3)],...
                  'CData',RGBim,'AlphaData', alpha,'ButtonDownFcn',{@SetPosition,'X',position(1),handles});
        end 
        % Zoom on IRM and draw lines
        volSize = size(stack(1).vol);
        corners = zeros(4,8);
        for k1 = 0:1
            for k2 = 0:1
                for k3 = 0:1
                    corners(:,k1*4+k2*2+k3+1) = stack(1).mat * [k1 * volSize(1) ; k2 * volSize(2) ; k3 * volSize(3) ; 1];
                end
            end
        end
        MIN = min(corners');
        MAX = max(corners');
      %  MIN = min(stack(1).mat * [0; 0; 0; 1],stack(1).mat * [volSize' ; 1]);
      %  MAX = max(stack(1).mat * [0; 0; 0; 1],stack(1).mat * [volSize' ; 1]);   
        axis(sp.z(j),'image',[MIN(1) MAX(1) MIN(2) MAX(2)]);
        axis(sp.y(j),'image',[MIN(1) MAX(1) MIN(3) MAX(3)]);    
        axis(sp.x(j),'image',[MIN(2) MAX(2) MIN(3) MAX(3)]);  
        set(sp.x(j),'Xdir','reverse','xticklabel',[],'yticklabel',[],'color','black');
        set(sp.y(j),'xticklabel',[],'yticklabel',[],'color','black');
        set(sp.z(j),'xticklabel',[],'yticklabel',[],'color','black');

        % Draw lines
        set(gcf,'currentAxes', sp.x(j));
        line([MIN(2) MAX(2)],[position(3) position(3)], 'hitTest','off');
        line([position(2) position(2)],[MIN(3) MAX(3)], 'hitTest','off');
        set(gcf,'currentAxes', sp.y(j));
        line([MIN(1) MAX(1)],[position(3) position(3)], 'hitTest','off');
        line([position(1) position(1)],[MIN(3) MAX(3)], 'hitTest','off');
        set(gcf,'currentAxes', sp.z(j));
        line([MIN(1) MAX(1)],[position(2) position(2)], 'hitTest','off');
        line([position(1) position(1)],[MIN(2) MAX(2)], 'hitTest','off');


    end
end

    function im2 = ApplyContrastBrightness(im)
       if stack(i).brightness < 0
            im = im*(1+stack(i).brightness);
        else
            im = im + (255-im)*stack(i).brightness;
        end
        im2 = uint16((int16(im)-128) * (tan((stack(i).contrast + 1) * pi/4)) + 128); 
    end

end



function SetPosition(src,event,type,offset,handles)
% src is image handle
% src.parent is the subplot
% src.parent.parent is the plotPanel
if strcmp(get(gcf,'selectionType'),'alt')
    f = figure('name','Full Screen');
    h = copyobj(get(src,'parent'),f);
    set(h,'position',[0 0 1 1]);
    images = get(h,'children');
    for i=1:length(get(h,'children'))
        set(images(i),'buttondownfcn',[]);
    end        
    delete(images(1));
    delete(images(2));
else
    pos = get(get(src,'parent'),'CurrentPoint');
    pos = pos(1,1:2);
    switch type
        case 'X'
            newPosition = [offset pos(1) pos(2)];
        case 'Y'
            newPosition = [pos(1) offset pos(2)];
        case 'Z'
            newPosition = [pos(1) pos(2) offset];
    end    
    mainFigureHandle = get(get(get(src,'parent'),'parent'),'parent');
    props.WindowButtonMotionFcn = get(mainFigureHandle,'WindowButtonMotionFcn');
    props.WindowButtonUpFcn = get(mainFigureHandle,'WindowButtonUpFcn');
    setappdata(mainFigureHandle,'TestGuiCallbacks',props);    
    set(mainFigureHandle,'WindowButtonMotionFcn',{@wbm,type,offset,get(src,'parent'),handles})
    set(mainFigureHandle,'WindowButtonUpFcn',{@wbu})
    SlicePlot(handles,newPosition);
end
end
% get the values and store them in the figure's appdata
% set the new values for the WindowButtonMotionFcn and
% WindowButtonUpFcn




% ---------------------------
function wbm(h,evd,type,offset,imageHandle,handles)
% executes while the mouse moves
 pos = get(imageHandle,'CurrentPoint');
    pos = pos(1,1:2);
    switch type
        case 'X'
            newPosition = [offset pos(1) pos(2)];
        case 'Y'
            newPosition = [pos(1) offset pos(2)];
        case 'Z'
            newPosition = [pos(1) pos(2) offset];
    end    
    SlicePlot(handles,newPosition);
end


% ---------------------------
function wbu(h,evd)
% executes when the mouse button is released
% get the properties and restore them
props = getappdata(h,'TestGuiCallbacks');
set(h,props);
end